Apparatus and method for sinogram restoration in computed tomography (CT) using adaptive filtering with deep learning (DL)

ABSTRACT

A method and apparatus is provided to reduce the noise in medical imaging by training a deep learning (DL) network to select the optimal parameters for a convolution kernel of an adaptive filter that is applied in the data domain. For example, in X-ray computed tomography (CT) the adaptive filter applies smoothing to a sinogram, and the optimal amount of the smoothing and orientation of the kernel (e.g., a bivariate Gaussian) can be determined on a pixel-by-pixel basis by applying a noisy sinogram to the DL network, which outputs the parameters of the filter (e.g., the orientation and variances of the Gaussian kernel). The DL network is trained using a training data set including target data (e.g., the gold standard) and input data. The input data can be sinograms generated by a low-dose CT scan, and the target data generated by a high-dose CT scan.

FIELD

This disclosure relates to reconstruction of medical images in whichdenoising of the images is performed using a deep-learning network basedon feature-aware training, and, more particularly, the denoising andartifact reduction can be performed on medical images including: (i)X-ray computed tomography (CT) images.

BACKGROUND

The background description provided herein is for the purpose ofgenerally presenting the context of the disclosure. Work of thepresently named inventors, to the extent the work is described in thisbackground section, as well as aspects of the description that cannototherwise qualify as prior art at the time of filing, are neitherexpressly nor impliedly admitted as prior art against the presentdisclosure.

Medical imaging produces images of the internal members of a patient'sbody. Examples of medical-imaging modalities include: X-ray radiography,X-ray computed tomography (CT), positron emission tomography (PET),single-photon emission CT (SPECT), fluoroscopy, and angiography. Oncethe images have been produced, a physician can use the images todiagnose a patient's injuries or diseases

X-ray CT systems and methods are widely used, particularly for medicalimaging and diagnosis. CT systems generally create images of one or moresectional slices through a subject's body. A radiation source, such asan X-ray source, irradiates the body from one side. At least onedetector on the opposite side of the body receives radiation transmittedthrough the body. The attenuation of the radiation that has passedthrough the body is measured by processing electrical signals receivedfrom the detector.

X-ray CT has found extensive clinical applications in cancer, heart, andbrain imaging. As CT has been increasingly used for a variety ofapplications including, e.g., cancer screening and pediatric imaging,there has arisen a push to reduce the radiation dose of clinical CTscans to become as low as reasonably achievable. For low-dose CT, theimage quality can be degraded by many factors, such as high quanta noisechallenge scanning geometry.

Although many cutting-edge technologies have been developed during toimprove low-dose CT image quality, better methods (e.g., faster, morerobust, and/or improved noise suppression) are desired to furthersuppress noise and generate clinical image quality with lower X-raydoses.

BRIEF DESCRIPTION OF THE DRAWINGS

A more complete understanding of this disclosure is provided byreference to the following detailed description when considered inconnection with the accompanying drawings, wherein:

FIG. 1 shows a flow diagram of a method of training a deep learning (DL)network and then using the DL network to denoise computed tomography(CT) data in the projection domain and then reconstruct an image fromthe denoised data, according to one implementation;

FIG. 2 shows a flow diagram of the process to denoise the projectiondomain data using the DL network, according to one implementation;

FIG. 3 shows a flow diagram of the process to train the DL network byiteratively adjusting coefficients of the DL network to optimize aloss-error function, according to one implementation;

FIG. 4A shows an example of a DL network, according to oneimplementation;

FIG. 4B shows an example of a type of DL network referred to as aconvolutional neural network (CNN), according to one implementation;

FIG. 5 shows diagram of a medical imaging system, according to oneimplementation;

FIG. 6 shows a schematic diagram of an X-ray CT scanner, according toone implementation;

FIG. 7A shows a perspective of a positron emission tomography (PET)scanner, according to one implementation; and

FIG. 7B shows a schematic diagram of the PET scanner, according to oneimplementation.

DETAILED DESCRIPTION

As discussed above, better methods for noise suppression are desirablefor many reasons, including that better noise suppression can enableclinical image quality at lower radiation doses. Additionally, methodsfor noise suppression can be improved by making them faster or morerobust. The methods described herein provide improved performance fornoise suppression by performing data-domain sinogram restoration usingdenoising filter that is specifically tailored to the given projectiondata. The denoising filter is specifically tailored to the givenprojection data using a deep learning network to determine for the givenprojection data the optimal parameters for a kernel of the denoisingfilter.

Accordingly, the methods described herein can advantageously reducecomputational time, hardware costs, and improve image quality formedical images, such as computed tomography (CT) images. Further, theexamples provided herein use CT imaging as an illustrative non-limitingexample. However, and the methods described herein can be used withother medical imaging modalities such as PET and SPECT, etc. by adaptingthe framework proposed herein. Accordingly, the discussion hereindiscloses and describes non-limiting examples of the present disclosure.As will be understood by those skilled in the art, the presentdisclosure may be embodied in other specific forms without departingfrom the spirit or essential characteristics thereof. Accordingly, thepresent disclosure is intended to be illustrative, but not limiting ofthe scope of the invention, as well as other claims. The disclosure,including any readily discernible variants of the teachings herein,defines, in part, the scope of the foregoing claim terminology such thatno inventive subject matter is dedicated to the public.

As discussed above, the methods described herein use a DL network (whichcan also be referred to as a neural network or an artificial neuralnetwork) that has been trained to select the parameters of a kernel of adenoising filter to optimize the sinogram restoration/denoising. Thus,the denoising filter is selectively tailored to denoise the particularprojection data. These kernel parameters selected by a DL network canresult in improved performance relative to related methods that also usedata-domain denoising techniques to improve the image quality of CTreconstructed images.

Instead of using a DL network, related methods use either adaptivefiltering or statistics-based mean estimation to select the kernelparameters. Compared to the methods described herein, these relatedmethods have their respective shortcomings/challenges. For example,statistical mean estimation methods require multiple iterations toobtain a converged solution, which can be computational intensive andtime consuming. Further, some of model parameters related to the meanestimation (e.g., parameter to control regularization strength) can onlybe selected manually, requiring user input.

Further, adaptive selection of the kernel parameters can be also bechallenging. For example, data-domain adaptive filtering method uses afilter having a small kernel (e.g., Gaussian function with a smallwidth/variance). The kernel is locally adapted to the data to smooth outnoise. The adaptive filters can be chosen from known kernel functions(e.g., triangle function and Gaussian function) with parameters such asthe variance of the kernel that are chosen empirically. Accordingly, theoptimal kernel design and parameter selection for adaptive filtering canbe difficult, especially in practical situations when the scanconditions are sometimes poor.

The methods described herein address the above-noted challenges withrelated methods by training a DL network to select the optimal kernelparameters for an adaptive filtering framework to perform sinogramrestoration/denoising in the data domain.

To illustrate, the methods described herein can use a smoothing filterto denoise projection/emission data in the data domain (as opposed todenoising in the image domain). A neural network can be trained todetermine the parameters for a kernel of the smoothing filter byapplying a sinogram to the neural network as an input. Further, thesmoothing filtering can be performed by convolving the sinogram with aGaussian kernel, the parameters of the Gaussian kernel (i.e., thewidth/variance of the Gaussian) that are determined by the neuralnetwork can be vary as a function of position/pixel within the sinogram.

Using a DL network to determine the kernel parameter for an adaptivesmoothing filter provides several advantages over the above-notedrelated methods.

First, the DL network can learn optimal strategies for adaptivelyfiltering in the data domain by training on actual data to minimize aloss function. This helps to reduce the effort of manual kernelselection which is often challenging in the related adaptive filteringmethod. Compared to ad-hoc or manually tuned adaptive filtering methods,the DL network can produce better performance (e.g., achieve a morefavorable tradeoff between noise reduction and resolution) by learningpatterns in the data that are too subtle or counterintuitive to discoverempirically without machine learning.

Second, the methods described herein are more flexible and robust thanthe related methods because the method for training the DL network isindependent of a particular type or source of noise.

Third, different loss function can be selected to achieve differentdenoising effects. For example, the loss function can be the

_(p)-norm of the difference between the target data and the result ofapplying the input data to the DL network. Different values of “p” inthe

_(p)-norm can be used to emphasize different aspects of the noise.Further, a weighting mask (e.g., based on the attenuation coefficient ofsignal intensity) can be applied on a pixel-by-pixel basis to thedifference between the target data and the result generated from theinput data. In certain implementations, rather than minimizing an

_(p)-norm of the difference between the target data and the result fromthe input data, the loss function can represent a similarity (e.g.,using a peak signal-to-noise ratio (PSNR) or a structural similarity(SSIM) index).

Fourth, the DL network can be fast computationally because only a smallnumber of parameters are required to define the kernel of the smoothingfunction.

Herein the phrase “data domain” is used to distinguish CT projectiondata from image data (i.e., reconstructed images) generated via CTreconstruction, which is in the “image domain.” That is, the “datadomain” is the projection data prior to reconstruction, and the “imagedomain” is the image data after reconstruction. The data domain can alsobe referred to as the “projection domain” and the “sinogram domain.” Ingeneral, the projection in the data domain has three dimensions, whichcan be referred to as view, channel, and segment, respectively. Thedimensions channel and segment correspond to the two directions/axes ofthe two-dimensional X-ray detector array, and the dimension viewcorresponds to the projection angle or projection view at which aprojection image is acquired. The denoising/smoothing filter can be athree-dimensional (3D) filter corresponding to all three dimensions ofthe data domain, but often the denoising/smoothing filter can be atwo-dimensional (2D) filter (e.g., a 2D convolution with a 2D Gaussiankernel). When the denoising/smoothing filter is a 2D filter, the filteris applied to only two of the three dimensions (e.g., to only the viewand channel dimensions or to only the channel and segment dimensions).

The term “data domain” is used instead of “sinogram domain” to avoid theerroneous interpretation that, because a sinogram plot often includesthe view dimension, the denoising/smoothing filter must be applied to aset of dimensions that includes the view dimension. Rather, the term“sinogram” is not limited to projection data that includes the viewdimension, but sinogram data can include, e.g., a projection image at a2D single view having only segment and channel dimensions. Accordingly,the data-domain denoising methods described herein are not limited toany particular combination of dimensions, but cover all possiblepermutations of dimensions for the denoising/smoothing filter (e.g., thedenoising/smoothing filter can be either a 2D or 3D filter applied toany combination of the view, channel, and segment dimensions). Forexample, the denoising/smoothing filter can being applied to (i) theview and channel dimensions, (ii) the channel and segment dimensions,(iii) the view and segment dimensions, and (iv) the view, channel, andsegment dimensions. Further, the methods described herein include theimplementation in which a first denoising/smoothing filter is applied toa first set of dimensions (e.g., the view and channel dimensions) togenerate once-filtered data, and then a second denoising/smoothingfilter is applied to the once-filtered data along a second set ofdimensions (e.g., the channel and segment dimensions) to generatetwice-filtered data.

X-ray CT is used as the primary illustrative example herein, but themethods described are also applicable to PET, SPECT, and fluoroscopy,for example. Although, the illustrative example describes applying asmoothing filter to sinograms from X-ray CT, the smoothing/denoisingcould also be applied to projection data from X-ray fluoroscopy acquiredat a single view angle because the smoothing/denoising is applied in thedata domain (i.e., not in the image domain). That is,smoothing/denoising can be applied in the data domain independent ofwhether the denoised projection data is later used for CTreconstruction.

Further, the filter applied in the data domain is not limited to being asmoothing filter. For example, the filter could include a convolutionkernel for edge enhancement or for artifact suppression. In certainimplementations, the convolution kernel and the parameters of theconvolution kernel can be selected to achieve one or more of denoising,edge enhancement, and artifact suppression. Nevertheless, the methodsdescribed herein are illustrated without loss of generality using thenon-limiting example of a smoothing/denoising filter to denoise theprojection data in the sinogram domain.

Referring now to the drawings, wherein like reference numerals designateidentical or corresponding parts throughout the several views, FIG. 1shows a flow diagram for a non-limiting example of a method that trainsand uses a DL neural network 170 to perform data-domain corrections ofCT projection data (e.g., sinogram restoration, denoising, smoothing,and/or artifact correction). Method 10, as illustrated in FIG. 1, usesDL network 170 to learn how to denoise raw data 105, and thenreconstructs a CT image from the denoised data. Method 10 includes twoparts: (i) an offline training process 150 and (ii) a medical imagingprocess 100. That is, process 150 trains the DL network 170, and process100 uses the trained DL network 170 to denoise projection data in thedata domain, thereby generating high-quality images 135 with reducednoise and artifacts.

In certain implementations, the network 170 includes a convolutionalneural network (CNN) in which series of convolution (conv), batchnormalization (BN), and rectified linear unit (ReLu) network layers areperformed.

The network 170 is trained using process 160. In process 160, a lossfunction is used to iteratively adjust/optimize parameters of the DLnetwork 170 (e.g., the parameters of the DL network 170 can includeweighting coefficients connecting network layers, and activationfunctions/potentials of nodes within the layers). The optimization ofthe network parameters continues until stopping criteria are satisfied(e.g., a stopping criterion can be whether the value of the lossfunction converged to a predefined threshold) to generate the trainednetwork 170.

The loss function compares target data 153 to an output acquired usingthe input data 157 and a current version of the DL network 170. Forexample, the input data can include noisy projection data in the datadomain for respective CT scans, and the target data can includelow-noise (i.e., less noise than the input data) projection data in thedata domain for the same respective CT scans. For example, each noisydataset for a given scan can form a pair with a correspondingnoise-reduced dataset for the given scan. In one implementation, thenoisy dataset is generated using a low-dose scan, and the noise-reduceddataset is generated using a high-dose scan. These scans can beperformed on a phantom, for example. Herein, to avoid the ambiguity ofusing the mass noun “data,” when referring to the projection data for agiven CT scan, the term “sinogram” is used herein to indicate the unitof projection data in the data domain corresponding to a complete CTscan.

Applying a noisy sinogram from the input data to the current version ofthe DL network 170 generates an output from the network that is supposedto be a denoised version of the noisy sinogram (i.e., a denoisedsinogram). The DL network 170 is trained by iteratively adjusting thenetwork coefficients in the DL network 170 to minimize the differencebetween the denoised sinogram output from the network 170 and thenoise-reduced sinogram from the target data 153. The training of thenetwork 170 is determined to be complete when the difference isminimized between network output and the target data, and whether or notthe difference has been sufficiently minimized is based on one or morepredetermined stopping criteria of process 160. Once the stoppingcriteria have been satisfied, the trained network 170 can then be storedand then later recalled to be used in the medical imaging process 100.

In method 10, a loss function is used to iteratively adjust networkcoefficients (e.g., weights and biases of convolutional and poolinglayers) of the DL network 170 until stopping criteria are satisfied(e.g., convergence of the parameters to a predefined threshold) togenerate the trained network 170. The loss function compareshigh-quality data 153 to results of a current version of the DL network170 to which input data 157 is applied.

As discussed above, CT image reconstruction is only one non-limitingillustrative example. Another example is positron emission tomography(PET) imaging. In the case PET imaging, sinograms can be generated forpositron emission data, and sinogram denoising/restoration can beapplied to the positron emission data using method 10. For example,method 10 includes training a DL network 170 and applying a low-quality(e.g., noisy) PET sinogram (i.e., raw data 105) to the trained network170 to generate a high-quality (e.g., denoised) PET sinogram.

In the case of PET imaging, the high- and low-quality data areaccumulated over scans having long and short time durations,respectively. In general, the signal-to-noise ratio (SNR) is smaller forsinograms accumulated over shorter time durations. Accordingly, thetarget data 153 (e.g., high-quality sinograms) can be generated usingall of the coincidence counts from a full-length PET scan to generatethe highest possible SNR for the sinogram. On the other hand, thelow-quality input data 157 can be generated using a partial subset ofcoincidence counts selected from the full dataset (e.g., using the datafrom only half of the full-length PET scan), resulting in a noisiersinogram (e.g., a √{square root over (2)} smaller SNR).

Returning to FIG. 1, process 100 is performed by obtaining raw data,e.g., by performing a CT scan to generate CT projections at a series ofview angles (i.e., a noisy sinogram). For example, the sinogram can beperformed using a low-dose CT scan to generate the raw data 105.

In step 110 of process 100, the raw data 105 is denoised by applying theraw data 105 to the trained DL network 170. The DL network 170 thenoutputs a denoised sinogram.

In step 120 of process 100, a CT image is reconstructed from thedenoised sinogram. Various methods can be used to reconstruct CT imagesfrom projection data, including filtered back-projection (FBP) andstatistical iterative reconstruction (IR) algorithms. In addition toFBP, other analytical methods can be used such as the Feldkamp DavisKress (FDK) method Adaptive Iterative Dose Reduction 3D (AIDR 3D)method. Compared to FBP reconstruction methods, IR methods can provideimproved image quality at reduced radiation doses.

One IR method performs unconstrained (or constrained) optimization tofind the argument p that minimizes the expression

${\underset{p}{argmin}\left\{ {{{{Ap} - \ell}}_{W}^{2} + {\beta\;{U(p)}}} \right\}},$wherein

is the projection data representing the logarithm of the X-ray intensityof projection images taken at a series of projection angles and p is areconstructed image of the X-ray attenuation for voxels/volume pixels(or two-dimensional pixels in a two-dimensional reconstructed image) inan image space. For the system matrix A, each matrix value a_(ij) (ibeing a row index and j being a column index) represents an overlapbetween the volume corresponding to voxel p_(j) and the X-raytrajectories corresponding to projection value

_(i). The data-fidelity term ∥Ap−

|_(w) ² is minimized when the forward projection A of the reconstructedimage p provides a good approximation to all measured projection images

. Thus, the data fidelity term is directed to solving the system matrixequation Ap=

, which expresses the Radon transform (i.e., projections) of variousrays from a source through an object OBJ in the space represented by pto X-ray detectors generating the values of

(e.g., X-ray projections through the three-dimensional object OBJ onto atwo-dimensional projection image

).

The notation ∥g∥_(w) ² signifies a weighted inner product of the formg^(T)Wg, wherein W is the weight matrix (e.g., expressing a reliabilityof trustworthiness of the projection data based on a pixel-by-pixelsignal-to-noise ratio). In other implementations, the weight matrix Wcan be replaced by an identity matrix. When the weight matrix W is usedin the data fidelity term, the above IR method is referred to as apenalized weighted least squares (PLWS) approach.

The function U(p) is a regularization term, and this term is directed atimposing one or more constraints (e.g., a total variation (TV)minimization constraint) which often have the effect of smoothing ordenoising the reconstructed image. The value β is a regularizationparameter is a value that weights the relative contributions of the datafidelity term and the regularization term.

In step 130 of process 100, additional image-domain denoising isperformed. This step is optional, and can be omitted in someimplementations.

Example denoising methods include linear smoothing filters, anisotropicdiffusion, non-local means, or nonlinear filters. Linear smoothingfilters remove noise by convolving the original image with a convolutionkernel that represents a low-pass filter or smoothing operation. Forexample, a Gaussian convolution kernel comprises elements determined bya Gaussian function. This convolution brings the value of each pixelinto closer agreement with the values of its neighbors. Anisotropicdiffusion removes noise while preserving sharp edges by evolving animage under a smoothing partial differential equation similar to theheat equation. A median filter is an example of a nonlinear filter and,if properly designed, a nonlinear filter can also preserve edges andavoid blurring. The median filter is one example of a rank-conditionedrank-selection (RCRS) filter, which can be applied to remove salt andpepper noise from an image without introducing significant blurringartifacts. Additionally, a filter using a total-variation (TV)minimization regularization term can be applied if imaged regionsupports an assumption of uniformity over large areas that are demarkedby sharp boundaries between the uniform areas. A TV filter is anotherexample of a nonlinear filter. Moreover, non-local means filtering is anexemplary method of determining denoised pixels using a weighted averageover similar patches within the images.

Finally, the reconstructed image 135 is generated, and the reconstructedimage 135 can be displayed to a used or stored for later use.

FIG. 2 show a flow diagram for a non-limiting example of step 110. Thesmoothing is performed in step 114, which uses a smoothing filter toremove noise by convolving the raw data with a kernel of a low-passfilter (e.g., a Gaussian). The shape of the kernel is not limited to aGaussian shape, but can be any shape (triangle, square, Blackman-Harris,Dolph-Chebyshev, Hann, Hamming, or other window shape). Further, thekernel does not have to be isotropic, but can be anisotropic. Further,the kernel can by dynamic as a function of position. That is, the kernelcan adapt to the local features of the raw data 105. The DL network 170is used to determine the filter parameters 113 that define the shape ofthe smoothing filter applied in step 114.

For example, if the smoothing filter is a bivariate Gaussian, then thefilter can be defined by three values: a first and second variance and arotation angle. Further, if the smoothing filter adapts as a function ofposition, then these three values can vary as a function of positionwithin the raw data 105.

In another example, if the smoothing filter is an isotropic Gaussian,then the filter can be defined by one value: the variance of theGaussian. And if the smoothing filter adapts as a function of position,then the variance of the Gaussian can vary as a function of positionwithin the raw data 105.

By applying the raw data 105 to the network 170, optimal filterparameters 113 can be determined as a function of position within theraw data 105. For example, a narrower kernel can be desirable in regionsof greater signal to preserve resolution and the noise is alreadyreduced in these regions. Further, sinograms can exhibit long skinnyregions in certain applications, and it might be desirable to havegreater smoothing along the long direction of a ridge than along thenarrow direction of the ridge. These are the types of characteristicsthat the DL network 170 through being trained in process 100.

In general, the raw data 105 can be either pre-log data (i.e.,proportional to the intensity of the X-ray radiation) or post-log data(i.e., proportional to the X-ray attenuation coefficient by taking thelogarithm of the X-ray intensity). Applying the raw data 105 to the DLnetwork 170 generate filter parameters, making the smoothing filteradaptive to the particular content of the raw data 105 (e.g., aspatially-varying rotation angle and variances). Then in step 114, theparameters generated by the DL network 170 are used to perform adaptivefiltering on the raw data 105, generating a restored/de-noised versionof the raw data 105.

In certain implementations, the DL network 170 used in step 112 is aconvolution neural network (CNN). The CNN can be a network that directlygenerates local small sized filters, e.g.,

$y_{i} = {\sum\limits_{j \in {{Neighbor}\mspace{14mu}{of}\mspace{14mu} i}}{w_{ij}x_{j}}}$wherein with w_(ij) is the filter on the ith pixel.

In certain implementations, the training data includes input data 157acquired via a low-dose scan and target data 153 acquired via ahigh-dose scan (i.e., a high-dose scan being any scan that uses agreater dose than the low-dose scan). Then the raw data 105 is acquiredusing a low-dose scan similar to that used to generate the input data157.

In certain implementations, the DL network 170 is a network thatgenerates kernel parameters for a Gaussian kernel; the kernel parameterdefining the variances and orientation angle based on the training usingthe training data. That is, a parametric filter is used (e.g., theGaussian kernel), and filtering is performed according to the parametersdetermined by applying the raw data 105 to the DL network 170.

In certain implementations, the target data 153 is higher quality rawdata or a sinogram (e.g., from a high dose scan).

Now a more detailed description of training a DL network is provided(e.g., process 160). This description is illustrated using the exampleof the target data 153 being the noise-reduced sinograms and the inputdata 157 being noisy sinograms.

FIG. 3 shows a flow diagram of one implementation of the trainingprocess 160. In process 160, input data 157 (e.g., noisy sinograms) andtarget data 153 (e.g., noise-reduced sinograms) are used as trainingdata to train a DL network 170, resulting in the DL network 170 beingoutput from step 319 of process 160. The offline DL training process 160trains the DL network 170 using a large number of input sinograms 157that are paired with corresponding target sinograms 153 to train the DLnetwork 170 to produce denoised sinograms resembling the targetsinograms 153 from the input sinograms 157.

In process 160, a set of training data is obtained, and the network 170is iteratively updated to reduce the error (e.g., the value produced bya loss function). The DL network infers the mapping implied by thetraining data, and the cost function produces an error value related tothe mismatch between the target sinograms 153 and the result produced byapplying a current incarnation of the DL network 170 to the inputsinograms 157. For example, in certain implementations, the costfunction can use the mean-squared error to minimize the average squarederror. In the case of a of multilayer perceptrons (MLP) neural network,the backpropagation algorithm can be used for training the network byminimizing the mean-squared-error-based cost function using a(stochastic) gradient descent method.

In step 316 of process 160, an initial guess is generated for thecoefficients of the DL network 170. For example, the initial guess canbe based on a priori knowledge of the region being imaged or one or moreexemplary denoising methods, edge-detection methods, and/or blobdetection methods. Additionally, the initial guess can be based on oneof a LeCun initialization, an Xavier initialization, and a Kaiminginitialization.

Steps 316 through 319 of process 160 provide a non-limiting example ofan optimization method for training the DL network 170.

An error is calculated (e.g., using a loss function or a cost function)to represent a measure of the difference (e.g., a distance measure)between the target sinograms 153 (i.e., ground truth) and inputsinograms 157 after applying a current version of the network 170. Theerror can be calculated using any known cost function or distancemeasure between the image data, including those cost functions describedabove. Further, in certain implementations the error/loss function canbe calculated using one or more of a hinge loss and a cross-entropyloss.

In certain implementations, the network 170 is trained usingbackpropagation. Backpropagation can be used for training neuralnetworks and is used in conjunction with gradient descent optimizationmethods. During a forward pass, the algorithm computes the network'spredictions based on the current parameters Θ. These predictions arethen input into the loss function, by which they are compared to thecorresponding ground truth labels (i.e., the high-quality data 153).During the backward pass, the model computes the gradient of the lossfunction with respect to the current parameters, after which theparameters are updated by taking a step of size of a predefined size inthe direction of minimized loss (e.g., in accelerated methods, such thatthe Nesterov momentum method and various adaptive methods, the step sizecan be selected to more quickly converge to optimize the loss function).

In certain implementations, the image processing in steps 112 and 114are considered as being the DL network for backpropagation. However,only the weighting coefficients in the CNN implemented in step 112 areallowed to be changed. That is, the weighting coefficients in the CNNcan be adjusted to generate better filter parameters 113, but the onlychanges to step 114 result from changes in filter parameters 113, whichoriginate outside of step 114. Step 114 itself remains fixed and thereare no changes internal to step 114. In this sense the entirety of steps112 and 114 can be considered as the DL network, even though only theweighting coefficients in the CNN of step 112 are being adjusted by thetraining process.

The optimization method by which the backprojection is performed can useone or more of gradient descent, batch gradient descent, stochasticgradient descent, and mini-batch stochastic gradient descent. Theforward and backwards passes can be performed incrementally through therespective layers of the network. In the forward pass, the executionstarts by feeding the inputs through the first layer, thus creating theoutput activations for the subsequent layer. This process is repeateduntil the loss function at the last layer is reached. During thebackward pass, the last layer computes the gradients with respect to itsown learnable parameters (if any) and also with respect to its owninput, which serves as the upstream derivatives for the previous layer.This process is repeated until the input layer is reached.

Returning to FIG. 3, step 317 of process 160 determines a change in theerror as a function of the change in the network can be calculated(e.g., an error gradient), and this change in the error can be used toselect a direction and step size for a subsequent change to theweights/coefficients of the DL network 170. Calculating the gradient ofthe error in this manner is consistent with certain implementations of agradient descent optimization method. In certain other implementations,this step can be omitted and/or substituted with another step inaccordance with another optimization algorithm (e.g., a non-gradientdescent optimization algorithm like simulated annealing or a geneticalgorithm), as would be understood by one of ordinary skill in the art.

In step 317 of process 160, a new set of coefficients are determined forthe DL network 170. For example, the weights/coefficients can be updatedusing the changed calculated in step 317, as in a gradient descentoptimization method or an over-relaxation acceleration method.

In step 318 of process 160, a new error value is calculated using theupdated weights/coefficients of the DL network 170.

In step 319, predefined stopping criteria are used to determine whetherthe training of the network is complete. For example, the predefinedstopping criteria can evaluate whether the new error and/or the totalnumber of iterations performed exceed predefined values. For example,the stopping criteria can be satisfied if either the new error fallsbelow a predefined threshold or if a maximum number of iterations isreached. When the stopping criteria is not satisfied the trainingprocess performed in process 160 will continue back to the start of theiterative loop by returning and repeating step 317 using the new weightsand coefficients (the iterative loop includes steps 317, 318, and 319).When the stopping criteria are satisfied the training process performedin process 160 is completed.

FIGS. 5A and 5B show various examples of the inter-connections betweenlayers in the DL network 170. The DL network 170 can include fullyconnected, convolutional, and the pooling layer, all of which areexplained below. In certain preferred implementations of the DL network170, convolutional layers are placed close to the input layer, whereasfully connected layers, which perform the high-level reasoning, areplace further down the architecture towards the loss function. Poolinglayers can be inserted after convolutions and proved a reductionlowering the spatial extent of the filters, and thus the amount oflearnable parameters. Activation functions are also incorporated intovarious layers to introduce nonlinearity and enable the network to learncomplex predictive relationships. The activation function can be asaturating activation functions (e.g., a sigmoid or hyperbolic tangentactivation function) or rectified activation function (e.g., theRectified Linear Unit (ReLU) applied in the first and second examplesdiscussed above). The layers of the DL network 170 can also incorporatebatch normalization, as also exemplified in the first and secondexamples discussed above.

FIG. 4A shows an example of a general artificial neural network (ANN)having N inputs, K hidden layers, and three outputs. Each layer is madeup of nodes (also called neurons), and each node performs a weighted sumof the inputs and compares the result of the weighted sum to a thresholdto generate an output. ANN s make up a class of functions for which themembers of the class are obtained by varying thresholds, connectionweights, or specifics of the architecture such as the number of nodesand/or their connectivity. The nodes in an ANN can be referred to asneurons (or as neuronal nodes), and the neurons can haveinter-connections between the different layers of the ANN system. Thesynapses (i.e., the connections between neurons) store values called“weights” (also interchangeably referred to as “coefficients” or“weighting coefficients”) that manipulate the data in the calculations.The outputs of the ANN depend on three types of parameters: (i) theinterconnection pattern between the different layers of neurons, (ii)the learning process for updating the weights of the interconnections,and (iii) the activation function that converts a neuron's weightedinput to its output activation.

Mathematically, a neuron's network function m(x) is defined as acomposition of other functions n_(i)(x), which can further be defined asa composition of other functions. This can be conveniently representedas a network structure, with arrows depicting the dependencies betweenvariables, as shown in FIG. 4A. For example, the ANN can use a nonlinearweighted sum, wherein m(x)=K(Σ_(i)w_(i)n_(i)(x)), where K (commonlyreferred to as the activation function) is some predefined function,such as the hyperbolic tangent.

In FIG. 4A (and similarly in FIG. 4B), the neurons (i.e., nodes) aredepicted by circles around a threshold function. For the non-limitingexample shown in FIG. 4A, the inputs are depicted as circles around alinear function, and the arrows indicate directed connections betweenneurons. In certain implementations, the DL network 170 is a feedforwardnetwork.

FIG. 4B shows a non-limiting example in which the DL network 170 is aconvolutional neural network (CNN). CNNs are type of ANN that hasbeneficial properties for image processing, and, therefore, havespecially relevancy for the applications of image denoising. CNNs usefeed-forward ANNs in which the connectivity pattern between neurons canrepresent convolutions in image processing. For example, CNNs can beused for image-processing optimization by using multiple layers of smallneuron collections which process portions of the input image, calledreceptive fields. The outputs of these collections can then tiled sothat they overlap, to obtain a better representation of the originalimage. This processing pattern can be repeated over multiple layershaving alternating convolution and pooling layers.

Following after a convolutional layer, a CNN can include local and/orglobal pooling layers, which combine the outputs of neuron clusters inthe convolution layers. Additionally, in certain implementations, theCNN can also include various combinations of convolutional and fullyconnected layers, with pointwise nonlinearity applied at the end of orafter each layer.

FIG. 5 illustrates an example embodiment of a medical-imaging system 40.The medical-imaging system 40 includes at least one scanning device 430;one or more image-generation devices 410, each of which is aspecially-configured computing device (e.g., a specially-configureddesktop computer, a specially-configured laptop computer, aspecially-configured server); and a display device 420.

The scanning device 430 is configured to acquire scan data by scanning aregion (e.g., area, volume, or slice) of an object (e.g., a patient).The scanning modality may be, for example, computed tomography (CT),positron emission tomography (PET), and/or single photon emission CT(SPECT). The one or more image-generation devices 410 obtain scan datafrom the scanning device 430 and generate an image of the region of theobject based on the scan data. After the one or more image-generationdevices 410 generate the image, the one or more image-generation devices410 send the image to the display device 420, which displays the image.

FIG. 6 illustrates in implementation in which the medical-imaging system40 includes a CT scanner. As shown in FIG. 7, a radiography gantry 500is illustrated from a side view and further includes an X-ray tube 501,an annular frame 502, and a multi-row or two-dimensional-array-typeX-ray detector 503. The X-ray tube 501 and X-ray detector 503 arediametrically mounted across an object OBJ on the annular frame 502,which is rotatably supported around a rotation axis RA.

The multi-slice X-ray CT apparatus further includes a high voltagegenerator 509 that generates a tube voltage applied to the X-ray tube501 through a slip ring 508 so that the X-ray tube 501 generates X-rays.The X-rays are emitted towards the object OBJ, whose cross sectionalarea is represented by a circle. For example, the X-ray tube 501 havingan average X-ray energy during a first scan that is less than an averageX-ray energy during a second scan. Thus, two or more scans can beobtained corresponding to different X-ray energies. The X-ray detector503 is located at an opposite side from the X-ray tube 501 across theobject OBJ for detecting the emitted X-rays that have transmittedthrough the object OBJ. The X-ray detector 503 further includesindividual detector elements or units.

The CT apparatus further includes other devices for processing thedetected signals from X-ray detector 503. A data acquisition circuit ora Data Acquisition System (DAS) 504 converts a signal output from theX-ray detector 503 for each channel into a voltage signal, amplifies thesignal, and further converts the signal into a digital signal.

The above-described data is sent to a preprocessing circuitry 506, whichis housed in a console outside the radiography gantry 500 through anon-contact data transmitter 505. The preprocessing circuitry 506performs certain corrections, such as sensitivity correction on the rawdata. A storage 512 stores the resultant data, which is also calledprojection data at a stage immediately before reconstruction processing.The storage 512 is connected to a processing circuitry 510 through adata/control bus 511, together with a reconstruction device 514, inputinterface 515, and display 516. The processing circuitry 510 controls acurrent regulator 513 that limits the current to a level sufficient fordriving the CT system.

The detectors are rotated and/or fixed with respect to the patient amongvarious generations of the CT scanner systems. In one implementation,the X-ray tube 501 and the X-ray detector 503 are diametrically mountedon the annular frame 502 and are rotated around the object OBJ as theannular frame 502 is rotated about the rotation axis RA.

The storage 512 can store the measurement value representative of theirradiance of the X-rays at the X-ray detector unit 503. Further, thestorage 512 can store a dedicated program for executing method 10.

The reconstruction circuitry 514 can execute various steps of method 10.Further, reconstruction circuitry 514 can execute pre-reconstructionprocessing image processing such as volume rendering processing andimage difference processing as needed.

The pre-reconstruction processing of the projection data performed bythe preprocessing circuitry 506 can include correcting for detectorcalibrations, detector nonlinearities, and polar effects, for example.

Post-reconstruction processing performed by the reconstruction circuitry514 can include filtering and smoothing the image, volume renderingprocessing, and image difference processing as needed. The imagereconstruction process can implement various steps of method 10. Thereconstruction circuitry 514 can use the memory to store, e.g.,projection data, reconstructed images, calibration data and parameters,and computer programs.

The various circuitry (e.g., the reconstruction circuitry 514 andpreprocessing circuitry 506) can include a CPU (processing circuitry)that can be implemented as discrete logic gates, as an ApplicationSpecific Integrated Circuit (ASIC), a Field Programmable Gate Array(FPGA) or other Complex Programmable Logic Device (CPLD). An FPGA orCPLD implementation may be coded in VHDL, Verilog, or any other hardwaredescription language and the code may be stored in an electronic memorydirectly within the FPGA or CPLD, or as a separate electronic memory.Further, the storage 512 can be non-volatile, such as ROM, EPROM, EEPROMor FLASH memory. The storage 512 can also be volatile, such as static ordynamic RAM, and a processor, such as a microcontroller ormicroprocessor, can be provided to manage the electronic memory as wellas the interaction between the FPGA or CPLD and the memory.

In one implementation, the reconstructed images can be displayed on adisplay 516. The display 516 can be an LCD display, CRT display, plasmadisplay, OLED, LED or any other display known in the art.

FIGS. 9A and 9B illustrates in implementation in which themedical-imaging system 40 includes a PET scanner that can implement themethod 10. The PET scanner includes a number of gamma-ray detectors(GRDs) (e.g., GRD1, GRD2, through GRDN) that are each configured asrectangular detector modules.

Each GRD can include a two-dimensional array of individual detectorcrystals, which absorb gamma radiation and emit scintillation photons.The scintillation photons can be detected by a two-dimensional array ofphotomultiplier tubes (PMTs) or silicon photomultipliers (SiPMs). Alight guide can be disposed between the array of detector crystals andthe photodetectors.

Each photodetector (e.g., PMT or SiPM) can produce an analog signal thatindicates when scintillation events occur, and an energy of the gammaray producing the detection event. Moreover, the photons emitted fromone detector crystal can be detected by more than one photodetector,and, based on the analog signal produced at each photodetector, thedetector crystal corresponding to the detection event can be determinedusing Anger logic and crystal decoding, for example.

FIG. 9B shows a schematic view of a PET scanner having gamma-ray(gamma-ray) photon counting detectors (GRDs) arranged to detectgamma-rays emitted from an object OBJ. The GRDs can measure the timing,position, and energy corresponding to each gamma-ray detection. In oneimplementation, the gamma-ray detectors are arranged in a ring, as shownin FIGS. 9A and 8B. The detector crystals can be scintillator crystals,which have individual scintillator elements arranged in atwo-dimensional array and the scintillator elements can be any knownscintillating material. The PMTs can be arranged such that light fromeach scintillator element is detected by multiple PMTs to enable Angerarithmetic and crystal decoding of scintillation event.

FIG. 9B shows an example of the arrangement of the PET scanner, in whichthe object OBJ to be imaged rests on a table 816 and the GRD modulesGRD1 through GRDN are arranged circumferentially around the object OBJand the table 816. The GRDs can be fixedly connected to a circularcomponent 820 that is fixedly connected to the gantry 840. The gantry840 houses many parts of the PET imager. The gantry 840 of the PETimager also includes an open aperture through which the object OBJ andthe table 816 can pass, and gamma-rays emitted in opposite directionsfrom the object OBJ due to an annihilation event can be detected by theGRDs and timing and energy information can be used to determinecoincidences for gamma-ray pairs.

In FIG. 9B, circuitry and hardware is also shown for acquiring, storing,processing, and distributing gamma-ray detection data. The circuitry andhardware include: a processor 870, a network controller 874, a memory878, and a data acquisition system (DAS) 876. The PET imager alsoincludes a data channel that routes detection measurement results fromthe GRDs to the DAS 876, a processor 870, a memory 878, and a networkcontroller 874. The data acquisition system 876 can control theacquisition, digitization, and routing of the detection data from thedetectors. In one implementation, the DAS 876 controls the movement ofthe bed 816. The processor 870 performs functions includingreconstructing images from the detection data, pre-reconstructionprocessing of the detection data, and post-reconstruction processing ofthe image data, as discussed herein.

The processor 870 can be configured to perform various steps of method10 described herein and variations thereof. The processor 870 caninclude a CPU that can be implemented as discrete logic gates, as anApplication Specific Integrated Circuit (ASIC), a Field ProgrammableGate Array (FPGA) or other Complex Programmable Logic Device (CPLD). AnFPGA or CPLD implementation may be coded in VHDL, Verilog, or any otherhardware description language and the code may be stored in anelectronic memory directly within the FPGA or CPLD, or as a separateelectronic memory. Further, the memory may be non-volatile, such as ROM,EPROM, EEPROM or FLASH memory. The memory can also be volatile, such asstatic or dynamic RAM, and a processor, such as a microcontroller ormicroprocessor, may be provided to manage the electronic memory as wellas the interaction between the FPGA or CPLD and the memory.

Alternatively, the CPU in the processor 870 can execute a computerprogram including a set of computer-readable instructions that performvarious steps of method 10, the program being stored in any of theabove-described non-transitory electronic memories and/or a hard diskdrive, CD, DVD, FLASH drive or any other known storage media. Further,the computer-readable instructions may be provided as a utilityapplication, background daemon, or component of an operating system, orcombination thereof, executing in conjunction with a processor, such asa Xeon processor from Intel of America or an Opteron processor from AMIDof America and an operating system, such as Microsoft VISTA, UNIX,Solaris, LINUX, Apple, MAC-OS and other operating systems known to thoseskilled in the art. Further, CPU can be implemented as multipleprocessors cooperatively working in parallel to perform theinstructions.

The memory 878 can be a hard disk drive, CD-ROM drive, DVD drive, FLASHdrive, RAM, ROM or any other electronic storage known in the art.

The network controller 874, such as an Intel Ethernet PRO networkinterface card from Intel Corporation of America, can interface betweenthe various parts of the PET imager. Additionally, the networkcontroller 874 can also interface with an external network. As can beappreciated, the external network can be a public network, such as theInternet, or a private network such as an LAN or WAN network, or anycombination thereof and can also include PSTN or ISDN sub-networks. Theexternal network can also be wired, such as an Ethernet network, or canbe wireless such as a cellular network including EDGE, 3G and 4Gwireless cellular systems. The wireless network can also be WiFi,Bluetooth, or any other wireless form of communication that is known.

While certain implementations have been described, these implementationshave been presented by way of example only, and are not intended tolimit the teachings of this disclosure. Indeed, the novel methods,apparatuses and systems described herein may be embodied in a variety ofother forms; furthermore, various omissions, substitutions and changesin the form of the methods, apparatuses and systems described herein maybe made without departing from the spirit of this disclosure.

The invention claimed is:
 1. An apparatus, comprising: circuitryconfigured to obtain radiation data representing an intensity ofradiation detected at a plurality of detectors, acquire a neural networkhaving weighting coefficients that have been trained to determineparameters of a filter, the parameters defining a shape function of aconvolution kernel of the filter, apply the radiation data to the neuralnetwork, and in response outputting parameters of the filter based onthe radiation data, filter the radiation data by applying the filterdefined by the output parameters from the neural network to generatefiltered radiation data, and reconstruct an image from the filteredradiation data.
 2. The apparatus according to claim 1, wherein thecircuitry is further configured to obtain the radiation data, whereinthe radiation data is a sinogram acquired by generating X-ray projectionimages at a series of projection views in which an angle between anX-ray source and the plurality of detectors is rotated relative to animage subject, and denoising the sinogram by applying the filter definedby the output parameters from the neural network generates a denoisedsinogram, wherein the filter is applied to a three-dimensional volume ofthe sinogram or to two-dimensional slices of the sinogram, which includea segment dimension and a channel dimension, a view dimension and thechannel dimension, and/or the segment dimension and a view dimension. 3.The apparatus according to claim 1, wherein the circuitry is furtherconfigured to obtain the radiation data, wherein the radiation data isone of X-ray computed tomography (CT) data, X-ray fluoroscopy data,gamma-ray positron emission tomography (PET), and single-photon emissionCT data (SPECT).
 4. The apparatus according to claim 2, wherein thecircuitry is further configured to reconstruct a computed tomography(CT) image from the denoised sinogram using an analytical reconstructionmethod.
 5. The apparatus according to claim 1, wherein the circuitry isfurther configured to filter the radiation data, wherein the filter is asmoothing filter, and filtering the radiation data denoises theradiation data.
 6. The apparatus according to claim 5, wherein thecircuitry is further configured to denoise the radiation data byapplying the smoothing filter, wherein the shape function of thesmoothing filter is a multi-variate, low-pass-filter kernel and theparameters of the smoothing filter comprise values defining widths andan orientation of the multi-variate, low-pass-filter kernel, and thevalues comprising the parameters vary as a function of pixel positionwithin the radiation data.
 7. The apparatus according to claim 6,wherein the circuitry is further configured to denoise the radiationdata by applying the smoothing filter, wherein the multi-variate,low-pass-filter kernel is a multi-variate Gaussian.
 8. The apparatusaccording to claim 1, wherein the circuitry is further configured tofilter the radiation data by applying the filter, wherein the filter isa two-dimensional filter applied to two-dimensional slices of theradiation data along a segment dimension and a channel dimension, a viewdimension and the channel dimension, or the view dimension and thesegment dimension.
 9. The apparatus according to claim 2, wherein thecircuitry is further configured to acquire the neural network, whereinthe weighting coefficients have been trained using training data thatincludes input data and target data, the input data comprising firsttraining sinograms acquired using a first radiation dose, and the targetdata comprising second training sinograms using a second radiation dosethat is greater than the first radiation dose.
 10. The apparatusaccording to claim 2, wherein the circuitry is further configured toobtain a training dataset comprising input sinograms paired withrespective target sinograms, the input sinograms exhibiting greaternoise than the corresponding target sinograms in the respective pairs,use the neural network to filter one of the input sinograms by applyingthe one of the input sinograms to the neural network and in responseoutputting parameters of the filter, and denoising the one of the inputsinograms by applying, to the one of the input sinograms, the filterdefined by the parameters output from the neural network, and therebygenerating a filtered sinogram, and train the neural network iterativelyadjusting the weighting coefficients to minimize a value of a lossfunction, the loss function measuring a disagreement between therespective target sinograms and the filtered sinograms corresponding tothe input sinograms.
 11. The apparatus according to claim 10, whereinthe circuitry is further configured to train the neural network whereinthe loss function includes a peak signal to noise ratio, a structuralsimilarity index, and/or an 1 p-norm of a difference between therespective target sinograms and the filtered sinograms corresponding tothe input sinograms.
 12. The apparatus according to claim 10, whereinthe circuitry is further configured to filter the one of the inputsinograms using the filter, wherein the filter is an adaptive filter inwhich the convolution kernel varies as a function of position within theone of the input sinograms thereby adapting to features represented inthe one of the input sinograms.
 13. A method, comprising: obtainingradiation data representing an intensity of radiation detected at aplurality of detectors; acquiring a neural network having weightingcoefficients that have been trained to determine parameters of a filterin response to an input, the parameters defining a shape function of aconvolution kernel of the filter; applying the radiation data to theneural network, and in response outputting the parameters of the filterbased on the radiation data; denoising the radiation data by applyingthe filter defined by the output parameters from the neural network togenerate filtered radiation data, and reconstructing an image from thefiltered radiation data.
 14. The method according to claim 13, whereinobtaining the radiation data includes that the radiation data is asinogram acquired by generating X-ray projection images at a series ofprojection views in which an angle between an X-ray source and theplurality of detectors is rotated relative to an image subject, anddenoising the sinogram by applying the filter defined by the outputparameters from the neural network generates a filtered sinogram,wherein the filter is applied to a three-dimensional volume of thesinogram or to two-dimensional slices of the sinogram, which include asegment dimension and a channel dimension, a view dimension and thechannel dimension, and/or the segment dimension and a view dimension.15. The method according to claim 13, wherein obtaining the radiationdata includes that the radiation data is one of X-ray computedtomography (CT) data, X-ray fluoroscopy data, gamma-ray positronemission tomography (PET), and single-photon emission CT data (SPECT).16. The method according to claim 13, wherein denoising the radiationdata by applying the filter further includes that the shape function ofthe filter is a multi-variate, low-pass-filter kernel and the parametersof the filter comprise values defining widths and an orientation of themulti-variate, low-pass-filter kernel, and the values comprising theparameters vary as a function of pixel position within the radiationdata.
 17. The method according to claim 14, wherein acquiring the neuralnetwork further includes that the weighting coefficients have beentrained using training data including input data and target data, theinput data comprising first training sinograms acquired using a firstradiation dose, and the target data comprising second training sinogramsacquired using a second radiation dose that is greater than the firstradiation dose.
 18. The method according to claim 14, further comprisingobtaining a training dataset comprising input sinograms paired withrespective target sinograms, the input sinograms exhibiting greaternoise than the corresponding target sinograms in the respective pairs,using the neural network to filter one of the input sinograms byapplying the one of the input sinograms to the neural network and inresponse outputting parameters of the filter, and denoising the one ofthe input sinograms by applying, to the one of the input sinograms, thefilter defined by the parameters output from the neural network, andthereby generating a denoisfiltered sinogram, and training the neuralnetwork iteratively adjusting the weighting coefficients to minimize avalue of a loss function, the loss function measuring a disagreementbetween the respective target sinograms and the denoised sinogramscorresponding to the input sinograms.
 19. The method according to claim18, wherein denoising the one of the input sinograms using the filterfurther includes that the filter is an adaptive filter in which theconvolution kernel varies as a function of position within the one ofthe input sinograms thereby adapting to features represented in the oneof the input sinograms.
 20. A non-transitory computer-readable storagemedium including executable instructions, which when executed bycircuitry, cause the circuitry to perform the method according to claim13.
 21. The apparatus according to claim 1, wherein the circuitry isconfigured to apply the radiation data to the neural network, and inresponse outputting the parameters changed for each position of theradiation data.